【直播】我的基因组26:GATK的realign对最后的结果影响到底有多大
GATK的realign步骤很耗费时间,使用起来也很繁琐,尤其是对初学者来说,耗费了大量的时间在学习,如何安装GATK软件,如何下载这个那个数据,但是这一步骤,但是不要只看别人这么做就盲目的跟着做,更要知其所以然。
realign对最后的variation的寻找真的意义重大吗?不见得,请看下面对结果具体的分析讲解。
因为我们这里还没有详细讲解vcf文件的格式,就先用现成的软件吧,snpEFF软件套装里面的SnpSift工具,具体安装教程见前面第5讲。
我也是第一次用,看软件说明书,需要我们把压缩的vcf先解压:
比较的第一个vcf文件是没有用GATK realign的,第二个是用了GATK realign的。但是因为我对它们的sample标记是一样的,所以在运行这个软件之前,我进入了用GATK realign的VCF文件里面手动编辑了sample标记。而且更改了文件名:
mv P_jmzeng.rmdup.realgn.bcftools.vcf realign.vcf
运行软件的命令是:
java -jar ~/biosoft/SnpEff/snpEff/SnpSift.jar concordance -v P_jmzeng.rmdup.bcftools.vcf realign.vcf 1>concordance.txt 2>SnpSift_Concordance.log
可能很多人不知道1> 2> 分别代表什么意思
在 shell 程式中,最常使用的 FD (file descriptor) 大概有三个, 分别是:
0 是一个文件描述符,表示标准输入(stdin)
1 是一个文件描述符,表示标准输出(stdout)
2 是一个文件描述符,表示标准错误(stderr)
这个软件运行很快,几分钟就好了!
首先看summary.txt文件,其实可以发现两个vcf文件差异并不大!
查看软件运行的日志(SnpSift_Concordance.log)可以知道,具体是哪些不match的
可以看到通过对GATK realign的的bam进行call snp的时候,多了很多indel的variation,而不是之前的snp情况了。而且还造成了一下alt的改变,
当然,也可以自己写脚本,来把两个vcf文件读进去,根据坐标来保存ref和alt,这样就能比较是否一致,我相信这个软件的原理也是一样的。
最后的结论,就是,GATK的realign步骤,对某些区域的variation情况的查询更精确,但是它改进的位点数量所占比例很有限,毕竟我这个vcf文件高达486万的variation,这个步骤仅仅是影响了三千多个位点的情况。
但是很多时候尤其是在临床应用,精确性是非常重要的,一个突变就可能意味着生与死的距离,但是科研不同,科研具有错误的机会。
接下来我们将重点讲解vcf格式文件里面所包含的信息。